Rate-Controlled  Constrained-Equilibrium  (RCCE)  Modeling  of 
Propulsive  Performance  of  Energetic  Materials  in  a  Hypersonic  Nozzle 


Final  Report 

Army  Grant  Number:  W911NF-08- 1-0444 


Mohamad  (Hameed)  Metghalchi 
Department  of  Mechanical  and  Industrial  Engineering 
Northeastern  University 
Boston,  MA 


October  2009 


Report  Documentation  Page 

Form  Approved 

OMB  No.  0704-0188 

Public  reporting  burden  for  the  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and 
maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information, 
including  suggestions  for  reducing  this  burden,  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington 

VA  22202-4302.  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  a  penalty  for  failing  to  comply  with  a  collection  of  information  if  it 
does  not  display  a  currently  valid  OMB  control  number. 

1.  REPORT  DATE 

2009  ^  REPORT  TYPE 

3.  DATES  COVERED 

00-00-2009  to  00-00-2009 

4.  TITLE  AND  SUBTITLE 

Rate- Controlled  Constrained-Equilibrium  (RCCE)  Modeling  of 
Propulsive  Performance  of  Energetic  Materials  in  a  Hypersonic  Nozzle 

5a.  CONTRACT  NUMBER 

W911NF-08-1-0444 

5b.  GRANT  NUMBER 

5c.  PROGRAM  ELEMENT  NUMBER 

6.  AUTHOR(S) 

5d.  PROJECT  NUMBER 

5e.  TASK  NUMBER 

5f.  WORK  UNIT  NUMBER 

7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS (ES) 

Department  of  Mechanical  and  Industrial  Engineering, Northeastern 
University, Boston, M  A, 021 15-5000 

8.  PERFORMING  ORGANIZATION  REPORT 

NUMBER 

;  55106-EG-II.l 

9.  SPONSORING/MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 

U.S.  Army  Research  Office,  P.O.  Box  12211,  Research  Triangle  Park, 
NC,  27709-2211 

10.  SPONSOR/MONITOR’S  ACRONYM(S) 

11.  SPONSOR/MONITOR’S  REPORT 
NUMBER(S) 

55106-EG-II.l 

12.  DISTRIBUTION/AVAILABILITY  STATEMENT 

Approved  for  public  release;  distribution  unlimited 

13.  SUPPLEMENTARY  NOTES 

14.  ABSTRACT 

15.  SUBJECT  TERMS 

16.  SECURITY  CLASSIFICATION  OF:  17.  LIMITATION 

_ _ _  OF  ABSTRACT 

18.  NUMBER  19a.  NAME  OF 

OF  PAGES  RESPONSIBLE  PERSON 

a.  REPORT  b.  ABSTRACT  c.  THIS  PAGE  Same  aS 

unclassified  unclassified  unclassified  Report  (SAR) 

20 

Standard  Form  298  (Rev.  8-98) 

Prescribed  by  ANSI  Std  Z39-18 


Abstract 


The  rate-controlled  constrained-equilibrium  (RCCE)  method  for  simplifying  the  kinetics  of 
complex  reacting  system  is  reviewed.  This  method  is  based  on  the  maximum  entropy  principle  of 
thermodynamics  and  involves  the  assumption  that  evolution  of  a  system  can  be  described  using 
few  set  of  slowly  changing  constraints  imposed  by  external  and  internal  dynamics  of  the  system. 
As  a  result,  the  number  of  differential  and  algebraic  equations  required  to  determine  the 
constrained-equilibrium  state  of  a  system  can  be  very  much  smaller  than  the  number  of  species. 
The  method  has  been  applied  to  the  oxidation  of  methane  in  a  constant  volume  adiabatic 
chamber.  The  results  of  using  12  constraints  (12  differential  equations)  are  in  perfect  agreement 
with  detailed  chemistry  using  29  species.  The  model,  then,  was  applied  to  the  expansion  of 
combustion  products  of  methane  during  the  power  stroke  of  an  internal  combustion  engine. 
Results  of  using  8  constraints  were  in  an  excellent  agreement  with  those  of  detailed  chemistry. 
Reacting  gas  mixture  through  a  nozzle  was  analyzed  using  detailed  chemistry,  and  properties 
such  as  temperature  and  pressure  as  a  function  of  location  were  determined.  The  combined 
analysis  of  flow  through  nozzle  and  rate-controlled  constrained-equilibrium  model  has  been 
developed  and  calculations  are  being  made  now. 


1.  Introduction 

The  development  of  models  for  describing  the  dynamic  evolution  of  chemically  reacting 
systems  is  a  fundamental  objective  of  chemical  kinetics.  The  conventional  approach  to  this 
problem  involves  first  specifying  the  state  and  species  variables  to  be  included  in  the  model, 
compiling  a  “full  set”  of  rate-equations  for  these  variables,  and  integrating  this  set  of  equations  to 
obtain  the  time-dependent  behavior  of  the  system.  Such  models  are  frequently  referred  to  as 
“Detailed  Kinetic  ModeF’s  (DKMs). 

The  problem  is  that  the  detailed  kinetics  of  C/H/O/N  molecules  can  easily  involve  hundreds 
of  chemical  species  and  isomers,  and  thousands  of  possible  reactions  even  for  system  containing 
only  Ci  molecules.  Clearly,  the  computational  effort  required  to  treat  such  systems  is  extremely 
large.  The  difficulties  are  compounded  when  considering  reacting  turbulent  flows,  where  the 
complexity  of  turbulence  is  added  to  that  of  the  chemistry. 

As  a  result  a  great  deal  of  effort  has  been  devoted  to  developing  methods  for  reducing  the 
size  of  DMKs.  Among  the  most  prominent  are:  Quasi  Steady  State  Approximation  (QSSA)  [1], 
Partial  Equilibrium  Approximation  [2],  Intrinsic  Low  Dimensional  Manifolds  (ILDM)  [3], 
Computational  Singular  Perturbation  (CSP)  [4],  Adaptive  Chemistry  [5],  Directed  Relation 
Graph  (DRG)  [6]  and  The  ICE-PIC  method  [7], 

An  alternative  approach,  originally  proposed  by  Keck  and  Gillespie  [8]  and  later  developed 
and  applied  by  Keck  and  co-workers  [9-11],  and  others  [12-15]  is  the  Rate-Controlled 
Constrained-Equilibrium  (RCCE)  method.  This  method  is  based  on  the  maximum  entropy 
principle  of  thermodynamics  and  involves  the  fundamental  assumption  that  slow  reactions  in  a 
complex  reacting  system  impose  constraints  on  its  composition  which  control  the  rate  at  which  it 
relaxes  to  chemical  equilibrium,  while  the  fast  reactions  equilibrate  the  system  subject  to  the 
constraints  imposed  by  the  slow  reactions.  Consequently,  the  system  relaxes  to  chemical 
equilibrium  through  a  sequence  of  constrained-equilibrium  states  at  a  rate  controlled  by  the 
slowly  changing  constraints. 

An  advantage  of  the  RCCE  method  is  that  it  does  not  require  a  large  DKM  as  a  starting  point. 
Instead,  one  starts  with  a  short  list  of  constraints  and  rate-limiting  reactions,  to  which  more  can 
be  added  to  improve  the  accuracy  to  the  desired  level.  If  the  only  constraints  are  those  imposed 
by  slowly  changing  state  variables,  the  RCCE  method  is  equivalent  to  a  local-chemical- 
equilibrium  (LTE)  calculation.  If  the  number  of  constraints  in  an  RCCE  model  is  identical  to  the 
number  of  species  in  a  DKM  model,  then  the  total  number  of  equations,  differential  plus 
algebraic,  to  be  solved  will  be  identical  and  the  results  will  be  similar  but  not  identical  to  those  of 
the  DKM  model.  The  DKM  will  include  only  the  specified  species  whereas  the  RCCE  model 
will  include  implicitly  all  species  which  can  be  made  from  the  elements  in  the  system. 


As  with  all  thermodynamic  systems,  the  number  of  constraints  necessary  to  describe  the 
dynamic  state  of  the  system  within  measurable  accuracy  can  be  very  much  smaller  than  the 


number  of  species  in  the  system.  Therefore,  fewer  equations  are  required  to  describe  the 
evolution  of  the  system.  A  further  advantage  is  that  only  rates  of  the  slow  constraint  changing 
reactions  are  needed  and  these  are  the  ones  most  likely  to  be  known.  Reactions  which  do  not 
change  any  constraint  are  in  equilibrium  and  need  not  be  specified.  It  should  be  emphasized  that 
the  successful  implementation  of  the  RCCE  method  depends  critically  on  the  constraints 
employed  and  a  knowledge  of  the  rates  of  the  constraint  changing  reactions  is  required.  When 
rates  are  uncertain  or  unknown,  the  Principle  of  Bayesian  Inference  suggests  that  the  best  choice 
is  to  assume  that  the  corresponding  reactions  are  in  constrained-equilibrium  (CE). 

The  purpose  of  this  project  was  to  investigate  the  application  of  RCCE  to  nozzle  flow 
analysis.  During  this  period,  evolution  of  a  reacting  system  in  a  constant  volume,  initially  far 
from  equilibrium  state,  was  modeled  successfully  using  RCCE  and  results  have  been  published 
in  Combustion  and  Flame  [16].  The  next  immediate  step  was  to  study  kinetic  effects  in  a  rapidly 
expanding  volume,  containing  a  mixture  initially  at  chemical  equilibrium.  The  model  considered 
was  expansion  of  combustion  products  in  a  cylinder  of  an  internal  combustion  engine.  The 
change  in  the  thermodynamic  state  due  to  expansion  causes  changes  in  chemical  composition. 
The  faster  the  expansion  the  lesser  the  time  to  adjust  to  the  new  thermodynamic  state  and  the 
more  pronounced  the  nonequilibrium  effects.  The  chemical  phenomena  observed  under  such 
conditions  are  expected  to  be  similar  to  those  in  a  rapidly  expanding  nozzle  flow.  The  results  of 
volume  expansion  have  been  published  in  International  Journal  of  Thermodynamics  [17].  The 
research  effort  has  been  then  followed  by  analyzing  flow  through  nozzle.  The  governing 
equations  using  detailed  kinetics  were  solved  and  results  are  presented.  The  combined  analysis  of 
flow  through  nozzle  and  rate-controlled  constrained-equilibrium  model  has  been  developed, 
shown  in  section  7.2,  and  calculations  are  being  made  now. 


2.  Rate-Controlled  Constrained-Equilibrium  (RCCE)  Method 

A  detailed  description  of  the  Rate-Controlled  Constrained-Equilibrium  (RCCE)  method  is 
given  in  reference  [18].  A  concise  summary  of  the  working  equations  for  chemically  reacting  gas 
mixtures  is  given  below. 

It  is  assumed  that  energy  exchange  reactions  are  sufficiently  fast  to  equilibrate  the 
translational,  rotational,  vibrational,  and  electronic  degrees  of  the  system  subject  to  constraints 
on  the  volume,  V,  and  the  total  energy  ,  E,  of  the  system.  Under  these  conditions,  the  energy  can 
be  written 

E=Et(T)N  (1) 

where  ET  (7)  is  the  transpose  of  the  species  molar  energy  vector  and  T  is  the  temperature  and  N 
is  the  vector  of  species  mole  numbers.  It  is  further  assumed  that  the  perfect  gas  model  applies 


and  the  species  constraints,  C,  can  be  expressed  as  a  linear  combination  of  the  species  mole 
numbers  in  the  form 


C=AN  (2) 

where,  A  is  the  nc  xnsp  constraint  matrix,  nc  is  the  number  of  constraints  ,  nsp  is  the  number  of 

species  and  C  and  N  are  column  vectors  of  length  nc  and  nsp .  Maximizing  the  Entropy,  S(E,  V, 
C),  subject  to  the  constraints  (2)  using  the  method  of  undetermined  Lagrange  multipliers,  [18], 
we  obtain  the  constrained-equilibrium  composition  of  the  system 

N!=(M/p)sx p(-/-/)  (3) 

where  M  is  the  total  mole  number,  p=MRT/V  is  the  pressure,  jj°  =  (h°  —Ts0)/  RT  is  the 
dimensionless  standard  Gibbs  free  energy  for  nsp  species,  and 

ft  i  (4) 

is  the  dimensionless  constrained-equilibrium  Gibbs  free  energy  of  the  species.  In  Eq.  (4),  A7 is 
the  transpose  of  the  constraint  matrix,  and  ^  is  the  vector  of  dimensionless  constraint  potentials 
(Lagrange  multipliers)  conjugate  to  the  constraint  vector,  C. 

Given  the  values  of  the  constraints,  C,  and  energy,  E,  of  the  system,  substitution  of  Eq.  (3) 
into  Eqs.  (1)  and  (2)  gives  a  set  of  nc+ 1  transcendental  equations  which  can  be  solved  for  the 
temperature,  T(Ey,Q  y°) ,  and  the  constraint  potentials,  y(E,V,C,  y(),  using  generalized 
equilibrium  codes  such  as  GNASA  [19]  or  GSTANJAN  [19].  Finally  substituting  Eq.  (4)  into 
Eq.  (3)  gives  the  constrained-equilibrium  composition  of  the  system,  N°  (E,V,C,y). 

2.1.  Rate-equations  for  the  constraints 

It  is  assumed  that  changes  in  the  chemical  composition  of  the  system  are  the  results  of 
chemical  reactions  of  the  type 

Y  K  <r^Y+  —  (5) 

where  X  is  the  species  vector,  y~  and  y/+are  nr  x  tfp  matrices  of  stoichiometric  coefficients  of 

reactants  and  products  respectively,  and  nr  is  the  number  of  reactions.  The  corresponding  rate- 
equations  for  the  species  can  be  written 

N=Vyr  (6) 

where  v  =  V+  —  V  >  r  =  r+~r~,  and  r+  and  r  are  the  forward  and  reverse  reaction  rate  column 
vectors  of  length  nr 


Differentiating  Eqs.  (1)  and  (2)  with  respect  to  time  and  using  Eq.  (6),  we  obtain  equations 
for  the  energy 

e=t£n+etn 

and  constraints 

C-AN-VBr 

where 

B=Ay  (9) 

is  an  nc  x«r  matrix  giving  the  change  of  constraints  due  to  the  n,  elementary  chemical  reactions 
among  species  and  Cv  =  QEI QT  is  the  molar  specific  heat  vector  at  constant  volume.  It  follows 

from  Eq.  (9)  that  a  reaction  k  for  which  all  Bik  are  zero  will  be  in  constrained-equilibrium  and  that  a 
constraint  i  for  which  all  B ik  are  zero  will  be  conserved.  The  latter  is  assumed  to  be  the  case  for  the 
elements. 

Given  equations  for  the  state-variables,  V(t)  and  E(t),  and  initial  values  for  the  species,  N( 0), 
Eqs.  (7)  and  (8)  can  be  integrated  in  stepwise  fashion  to  obtain  the  temperature,  T(t)  and  species 
constraints,  C(t).  At  each  time  step,  a  generalized  equilibrium  code,  such  as  GNASA  or  GSTANJAN 
previously  cited,  must  be  used  to  determine  the  temperature,  T(E,V,Q  j/ ),  and  constrained- 

equilibrium  composition,  Et  ( E,  V,  C,  jj/’  ).  These,  in  turn,  can  be  used  to  evaluate  the  reaction  rates 

r(T,V,  jj°  ,  y  )  required  for  the  next  step.  Note  that  only  the  rates  of  reactions  which  change 

constraints,  i.e.  those  for  which  Bik  ^0,  are  required  for  RCCE  calculations.  All  other  reactions  are  in 
constrained-equilibrium  and  need  not  be  specified. 


(7) 

(8) 


3.  Rate  equations  for  the  constraint  potentials 

Although  direct  integration  of  the  rate-equations  for  the  constraints  is  relative  straight  forward 
and  simple  to  implement,  it  has  proved  to  be  relatively  inefficient  and  time  consuming  due  to  the 
slowness  of  the  constrained-equilibrium  codes  currently  available  [19].  An  alternative  method,  first 
proposed  by  Keck  [18]  and  implemented  in  later  works  [19,20]  with  colleagues,  is  the  direct 
integration  of  the  rate-equations  for  the  constraint-potentials.  This  method  has  also  recently  been 
investigated  by  Tang  and  Pope  [5]  and  Jones  and  Rigopoulos  [6]. 

Differentiating  Eq.  (3)  with  respect  to  time  and  substituting  the  result  into  Eq.  (7)  and  (8) 
yields  the  nc+l  implicit  equations  for  the  energy, 


(10) 


=r y  y 

where 

Dp  =  ^y,ijEjNj  (10a) 

Dv  =^EjNCj  (10b) 

Dr  =^\CV]T +£*  / RT)Ncj  (10c) 

and  the  constraint-potentials, 

£Tz-Qv^-QT^+ii= 0  (11) 

where 

Cnk=^yijakjNcj  (11a) 

(lib) 

Cn  =  ^OijEjNj  /RT  (lie) 


In  this  case,  given  equations  for  the  state  variables,  V(t)  and  E(t),  and  the  initial  temperature, 
T(0),  and  constraint-potentials,  y(0),  Eqs.  (10)  and  (11)  can  be  integrated  using  implicit  ODE 

integration  routines  such  as  DASSL  [7]  to  obtain  the  temperature,  T(t),  and  constraint  potentials, 
jXt) .  The  constrained-equilibrium  composition,  P/  (E,  V,  t),  of  the  system  can  then  be  determined 

using  Eq.  (3).  The  number  of  unknowns  is  reduced  from  the  nsp+l  in  a  DKM  calculation  to  nc+l 
in  the  RCCE  calculation.  As  previously  noted,  only  the  rate  constants  for  those  reactions  which 
change  constraints,  i.e.  Bik  ±  0,  are  needed.  Note  that,  once  the  constraint-potentials  have  been 
determined,  the  constrained-equilibrium  concentration  of  any  species  for  which  the  standard 
Gibbs  free  energy  is  known  can  be  calculated  whether  or  not  it  is  explicitly  included  in  the 
species  list.  It  also  follows  that  all  species  are  implicitly  included  in  the  RCCE  rate-equations. 


4.  Selection  of  constraints 

The  careful  selection  of  constraints  is  the  key  to  the  success  of  the  RCCE  method.  Among  the 
general  requirements  for  the  constraints  are  that  they  must  a)  be  linearly  independent 
combinations  of  the  species  mole  numbers,  b)  hold  the  system  in  the  specified  initial  state,  c) 
prevent  global  reactions  in  which  reactants  or  intermediates  go  directly  to  products  and,  d) 


determine  the  energy  and  entropy  of  the  system  within  experimental  accuracy.  In  addition,  they 
should  reflect  whatever  information  is  available  about  rate-limiting  reactions  which  control  the 
evolution  of  the  system  on  the  time  scale  of  interest. 

In  the  present  work,  the  focus  is  on  applications  of  the  RCCE  method  to  chemically  reacting 
gas  phase  mixtures.  In  the  temperature  and  pressure  range  of  interest,  the  rates  of  nuclear  and 
ionization  reactions  are  negligible  compared  to  those  for  chemical  reactions  and  the  fixed 
constraints  are  the  neutral  elements  of  hydrogen,  carbon,  oxygen,  nitrogen,...,  designated  by  EH, 
EC,  EO,EN,... 

Under  these  conditions,  the  slowest  reactions  controlling  the  chemical  composition  are  three- 
body  dissociation/recombination  reactions  and  reactions  which  make  and  break  valence  bonds. 
Such  reactions  are  slow  in  the  endothermic  direction  because  of  the  high  activation  energies 
required,  and  in  the  exothermic  direction  because  of  small  three-body  collision  rates  and  small 
radical  concentrations.  They  impose  slowly  varying  time-dependent  constraints  on  the  number  of 
moles,  M,  of  gas  and  the  free  valence,  FV,  of  the  system  respectively.  A  finite  value  of  FV  is  a 
necessary  condition  for  chain  branching  chemical  reactions  to  proceed. 

A  third  important  time-dependent  constraint,  imposed  by  slow  OO  bond-breaking  reactions,  is 
the  free  -  oxygen ,  FO,  defined  as  any  oxygen  atom  not  directly  bound  to  another  oxygen  atom.  An 
increase  in  FO  is  a  necessary  condition  for  the  formation  of  the  major  reaction  products  of 
hydrocarbon  oxidation,  H2O,  CO2  and  CO. 

Two  additional  time-dependent  constraints  which  slightly  improve  the  agreement  between 
RCCE  and  DKM  calculations  under  some  conditions  are:  OHO  =  OH+O  and  DCO= 

HCO+CO.  The  OHO  constraint  is  a  consequence  of  the  relatively  slow  constraint  changing 
reaction  RH+OH  .^lUO+R  coupled  with  the  fast  reaction  RH+0=OH+R  which  equilibrates  OH 
and  O.  The  DCO  constraint  is  a  consequence  of  the  slow  spin-forbidden  reaction  CO+HO2  <-> 
CO2+OH  coupled  with  the  fast  reaction  HC0+02=C0+H02  which  equilibrates  HCO  and  CO. 

For  systems  involving  the  elements  C,  H  and  O  the  8  constraints  EH,  EO,  EC,  M,  FV,  FO, 
OHO,  and  DCO  are  independent  of  the  initial  reactants  and  may,  therefore,  be  considered 
“universal”  constraints.  Along  with  the  equilibrium  reactions, 

H2+0=OH+H  (CE1) 

H2+H00=H202+H  (CE2) 

HC0+02=C0+H02  (CE3) 

they  are  sufficient  to  determine  the  constrained-equilibrium  mole  fractions  of  the  11  major 
hydrocarbon  combustion  products  H,  O,  OH,  HO2,  H2,  O2,  H2O,  H2O2,  HCO,  CO  and  CO2 
under  both  high  and  low  temperature  conditions. 

In  the  present  investigation  of  Ci -hydrocarbon  oxidation,  4  additional  fuel-dependent 
constraints  have  been  used.  The  first  is  a  constraint  on  the  fuel,  FU,  imposed  by  slow  hydrogen  - 
abstraction  reactions  of  the  type  FU+O2  <->FR+H02  and  even  slower  dissociation/recombination 
of  the  type  AB+M  ^A+B+M.  This  constraint  is  necessary  to  hold  the  system  in  its  initial  state. 


The  second  is  a  constraint  on  fuel  radicals,  FR,  which  is  necessary  to  prevent  the  equilibration  of 
forbidden  exothermic  global  reactions  such  as  CH3+202+2H20=C02+2H202+H2+H  which 
would  otherwise  convert  fuel  radicals  directly  to  major  products.  The  third  is  a  constraint  on 
alkylperoxide s,  APO  =  CH3OOH+CH3OO+CH2OOH,  imposed  by  slow  reactions  which  convert 
APO  to  hydroperoxides  coupled  with  fast  reactions  which  equilibrate  the  species  comprising 
APO,  and  the  fourth  is  a  constraint  on  alcohol  plus  formaldehyde,  ALCD  = 
CH3OH+CH3O+CH2OH+CH2O  imposed  by  relatively  slow  reactions  which  generate/remove 
ALCD  coupled  with  fast  reactions  which  equilibrate  the  species  comprising  ALCD. 


5.  CH4  Results  for  133  Reaction  Set 


Time-dependent  temperature  profiles  of  stoichiometric  mixtures  of  methane  and  oxygen  at 
initial  temperatures  of  900K  and  1500K  and  different  initial  pressures  are  shown  in  Fig.  1.  The 
12  constraints  used  are  listed  in  table  1.  Under  this  set  of  constraints,  only  102  of  the  reactions  in 
the  full  set  of  133  reactions  change  constraints  and  are  therefore  required.  The  remaining  31  are 
in  constrained-equilibrium  and  are  not  needed  in  RCCE  calculations. 


Log  t  (sec) 


Figure  1 :  Comparison  of  RCCE  and  DKM 
temperature  vs  time  profiles  for  a 
stoichiometric  mixture  of  CH4/O2  at 
different  initial  pressures  and  initial 
temperatures  of  900  K  and  1500  K.  (Time 
in  sec) 


Figure  2:  Comparison  of  RCCE  and  DKM 
temperature  vs  time  profiles  for  CH4/O2  mixtures 
with  different  equivalence  ratios  at  initial 
temperatures  and  pressures  of  900  K  and  100 
atm  (a)  and  1500  K  and  1  atm  (b).  (time  in  sec) 


It  can  be  seen  that  the  agreement  with  DKM  calculations  is  excellent  over  the  entire  range  of 
pressure  and  temperature  covered.  Note  that  the  temperature  overshoot  at  low  pressures  due  to 


slow  three-body  recombination  and  dissociation  reactions  is  well  reproduced  by  the  constraint  M 
on  the  total  moles. 

Additional  results  are  shown  on  log-log  plots  in  Fig.  2  for  the  initial  conditions  900  K  and 
100  atm,  where  the  dominant  radicals  are  H02~CH3  at  early  times  and  CH300  at  late  times  and 
in  Fig.  3  for  the  initial  conditions  1500  K  and  latm,  where  the  dominant  radicals  are 
H02~CH300  at  all  times.  As  can  be  seen  in  Figs  2a  and  3a,  the  temperature  first  decreases  due 
to  the  fact  that  the  initiation  reactions  are  all  endothermic  then  later  increases  as  exothermic 
reactions  become  important.  The  agreement  between  RCCE  and  DKM  calculations,  especially 
with  regard  to  the  time  at  which  the  temperature  difference  becomes  positive,  is  excellent. 

The  fixed  elemental  constraints  and  the  most  important  time-dependent  constraints  M,  FV, 
FO,  FU  and  FR  are  shown  in  Figs  2d  and  3d  and  the  mole  fractions  of  the  major  species  are 
shown  in  Figs  2e-h  and  3e-h.  It  can  be  seen  that  overall  agreement  is  very  good.  A  detailed 
presentation  of  this  research  has  been  published  in  Combustion  and  Flame  [16]. 


6.  Power  Stroke  of  and  Internal  Combustion  Engine 


The  physical  model  considered  for  this  part  of  the  research  is  shown  in  figure  4.  It  is  assumed 
that  the  power  stroke  begins  with  the  combustion  products  at  a  complete  chemical  equilibrium 
state  and  that  expansion  occurs  so  rapidly  that  heat  conduction  does  not  occur  through  the  piston 
wall.  It  is  further  assumed  that  the  gas  composition  is  homogeneous  at  any  time  during 
expansion.  Volume  is  assumed  to  be  a  prescribed  function  of  time. 
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Figure  4:  Schematic  of  the  physical  model 


Figure  3:  Comparison  of  various  RCCE  and  DKM  time  profiles  for  stoichiometric  CH4/02 
mixtures  at  initial  temperature  and  pressure  of  900  K  and  100  atm.  (time  in  sec) 
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Figure  4:  Comparison  of  various  RCCE  and  DKM  time  profiles  for  stoichiometric  CH4/02 
mixtures  at  initial  temperature  and  pressure  of  1500  K  and  1  atm.  (time  in  sec) 


Constraint 

Definition  of  the  Constraint 

1 

EC 

Elemental  carbon 

2 

EO 

Elemental  oxygen 

3 

EH 

Elemental  hydrogen 

4 

M 

Total  number  of  moles 

5 

FV 

Moles  of  free  valance  (any  unpaired  valence  electron) 

6 

FO 

Moles  of  free  oxygen  (any  oxygen  not  directly  attached  to  another  oxygen) 

7 

OHO 

Moles  of  water  radicals  (O+OH) 

8 

DCO 

Moles  of  HCO+CO 

9 

FU 

Moles  of  fuel  molecules 

10 

FR 

Moles  of  fuel  radicals 

11 

APO 

Moles  of  AlkyPerOxydes  (CH300+CH300H+CH200H) 

12 

ALCD 

Moles  of  Alcohols+Aldehydes  (CH30+CH30H+CH20H+CH20) 

The  general  formulation  remains  the  same.  However,  the  time  rate  of  change  of  volume  in 
equations  (10)  and  (11)  is  no  longer  zero  in  this  case.  The  constraints  used  in  this  study  include 
EH,  EO,  EC,  EN  (elemental  nitrogen),  M,  FV,  FO  and  DCO  from  table  1. 

According  to  the  Le  Chatelier  principle,  when  a  highly  dissociated  mixture  undergoes  an 
interaction  with  the  environment  which  lowers  the  gas  temperature  and  density,  the  internal 
dynamics  shift  in  the  exothermic  direction  to  minimize  the  cooling  effect  of  interaction.  As  a 
result,  three  body  recombination  reactions  become  an  important  part  of  the  energy  restoration 
process.  The  total  number  of  particles  in  the  system  does  not  change  unless  one  such  reaction 
occurs.  Figures  5  and  6  represent  comparison  of  the  detailed  kinetics  against  RCCE  when 
constraints  are  added  on  at  a  time.  RCCE  is  able  to  handle  cases  as  different  as  frozen 
equilibrium,  fully  rate  controlled  and  local  thermodynamic  equilibrium  (L.T.E)  with  the  same 
ease.  The  result  of  this  part  of  the  research  has  been  published  in  International  Journal  of 
Thermodynamics  [17]. 

7.  Governing  Equations  for  an  Expanding  Nozzle  Flow 
7.1.  Detailed  Kinetic  Formulation 

Assuming  an  inviscid,  adiabatic  and  steady  flow,  it  is  easy  to  show  that  the  enthalpy  is 
conserved  and  the  governing  equations  take  the  following  form: 

Overall  Continuity 

ftiAx)  —m—cte  (12) 


Frozen  Equilibrium 


Figure  5:  Temperature  profiles  under  different 
sets  of  constraints  compared  with  the  predictions 
of  detailed  kinetics  at  an  engine  speed  of  3000 
rpm 
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Figure  6:  Constraint-dependence  of  CO  profiles 
compared  with  the  predictions  of  detailed 
kinetics  at  an  engine  speed  of  3000  rpm 


Momentum  Equation 


Species  Mass  Balance 


Energy  Equation 


Equation  of  State 


du  dp  n 
U  dx+  dx~ 


dx 

dx 


+u 


du  q 
dx  ~ 


(13) 


(14) 


(15) 


(16) 


where  p,  u,  A(x) ,  p  and  T  are  density,  convective  velocity,  nozzle  area,  pressure  and 
temperature  respectively.  Also,  Yt,  1r\,  Cpi  and  A^s/) represent,  respectively,  the  mass  fraction, 

total  enthalpy  (formation  +  sensible),  heat  capacity  at  constant  pressure  of  species  i  and  number 
of  species.  Qi  is  the  net  molar  rate  of  consumption  of  species  i  and  M(  is  its  molecular  mass. 

The  mixture  molecular  mass  Mis  also  determined  as  follows: 


(17) 


M  ~ 


Note  that  the  energy  equation  (15)  can  be  further  simplified  by  replacing  for 


— —  from  equation 
ax 


(14)  as: 


The  initial  condition  is  taken  to  be  a  C/H/O  system  at  chemical  equilibrium.  Such  initial 
condition  can  be  easily  obtained  using  STANJAN  program.  The  flow  is  subsonic  up  to  the  throat 
section.  Therefore,  it  is  a  reasonable  to  assume  that  the  mixture  is  in  chemical  equilibrium. 
Departures  from  equilibrium  happen  due  to  expansion  in  the  divergent  section  as  a  result  of 
insufficient  residence  time.  In  the  current  study  we  consider  the  expansion  downstream  of  the 
throat  section  and  specify  an  equilibrium  mixture  at  this  section.  Equations  (13)  through  (18)  are 
then  solved  under  an  area  profile  presented  in  section  7.3.  The  initial  temperature  is  2000  K  and 
the  initial  pressure  is  20  atm  and  the  initial  velocity  is  900  m/sec.  Figs  7  and  8  show  the 
temperature  and  pressure  profiles  within  the  divergent  section  of  the  nozzle.  The  exit-to-throat 
area  ratio  is  assumed  10,  the  nozzle  length  is  assumed  2m,  and  the  throat  area  lm2. 


Figure7:Temperature  profiles  along  the  axis  of  Figure  8:  Pressure  profiles  along  the  axis  of  nozzle, 
nozzle. 


7.2.  RCCE  Formulation 

Equations  (12)  and  (13)  remain  invariant  under  RCCE  transformations.  In  what  follows  the 
RCCE  form  of  the  species  continuity,  the  energy  equation  and  the  equation  of  state  are  derived. 


Species  Mass  Balance 

Multiplying  equation  (14)  through  by  A(x)and  then  using  the  continuity  equation,  and  also  using 
the  relation  p  =[Ni]Mi ,  one  can  easily  obtain  the  following  form  for  the  species  equation: 


rAT1du  d[N: ]  VAT/1  dA\ 


(19) 


dx  dx 
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Adx  j 


The  RCCE  distribution  of  concentrations  is  also  given  as  [1]: 


Nr 


(20) 


(21) 


[A^]=QexpC-gflyp) 

Differentiating  equation  (20)  with  respect  to  X yields: 

(Emyi-™bdl 


V  V  , 

Now,  multiplying  equation  (19)  by  Clik ,  summing  over  species  index  i  and  then  using  the  equation 
(21),  we  will  have: 
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Energy  Equation: 

The  energy  equation  in  this  case  is  equivalent  with  the  conservation  of  stagnation  enthalpy,  that 
is: 


//0  =/t+^M2  =cte 


(23) 


where 


The  energy  equation  can  then  be  written  in  the  following  form: 


N  N 

1  ysp  1 ysp 


(25) 


Differentiating  equation  (25)  with  respect  to  X and  using  equation  (2)  results  in  the  following 
form  for  energy  equation: 


■YE id^+Er^f+fp^,+du2-Hji’=0 

4h  y  dx  dx  H  dx  2  ax 


(25) 


where 


and  H0  is  the  stagnation  enthalpy  per  unit  mass  at  nozzle  throat,  which  is  constant  everywhere. 


Equation  of  State: 

Re-writing  equation  (16)  using  the  relation  p  =[Nj]Mi in  equation  (17)  results  in  the  following 
familiar  form  for  the  equation  of  state: 


AT 


p=^[  njkt 

Differentiating  equation  (26)  with  respect  to  Xand  then  using  equation  (21)  yields: 


(26) 


where 


dT  dp  q 
dx  dx~ 


(27) 


s= 


-K 


7.3  Nozzle  Area  Profile 


Currently,  the  following  profile  has  been  used: 


A(x)  —ax:  +Ivc  +cx+d  (28) 

The  unknown  coefficients  a,  b,C,d  are  determined  such  that  the  following  geometrical  constraints 
are  satisfied: 


A{ xj)=Ah’ 

A^x=L—  ^end> 

—  I  -0 

dx  xJ)~  ’ 

—  I  -0 

The  nozzle  geometry,  under  these  boundary  conditions,  becomes: 

^ =-2(4  -if- +3(4  -l/ff  +1 

Ah 


Where 


(29) 


(30) 
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